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Abstract 

We use density-functional theory to study the structure of two-dimensional defects inside a 
circular nematic nanocavity. The density, nematic order parameter, and director fields, as well as 
the defect core energy and core radius, are obtained in a thermodynamically consistent way for 
defects with topological charge k = +1 (with radial and tangential symmetries) and k = +1/2. An 
independent calculation of the fluid elastic constants, within the same theory, allows us to connect 
with the local free-energy density predicted by elastic theory, which in turn provides a criterion 
to define a defect core boundary and a defect core free energy for the two types of defects. The 
radial and tangential defects turn out to have very different properties, a feature that a previous 
Maier-Saupe theory could not account for due to the simplified nature of the interactions -which 
caused all elastic constants to be equal. In the case with two k = +1/2 defects in the cavity, the 
elastic regime cannot be reached due to the small radii of the cavities considered, but some trends 
can already be obtained. 

PACS numbers: Valid PACS appear here 
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I. INTRODUCTION 



The analysis of defects in liquid crystals is very important from many points of view. In 
liquid-crystal applications, defects play a crucial role in governing display-cell operation. 
Also, there are interesting theoretical issues in different areas of physics concerning defects 
and the stabilisation of defects has been observed and analysed in computer simulations 
[10 , [5|. A defect is a singularity in the director field of the liquid crystal . Local 
properties of the liquid crystal, e.g. the nematic order parameter, asymptotically relax to 
values of the bulk material far from the singularity but, in its immediate neighbourhood, 
properties undergo abrupt (i.e. within molecular lengths) changes; this region somehow 
defines microscopically a boundary for the so-called defect core. 

Beyond the defect core, variations are smooth, so that the macroscopic elastic theory 
of Frank [8j can be used, together with some assumptions about defect core energies and 
radii. Very often core energies are simply ignored. It would be desirable to have estimations 
of these properties based on more microscopic approaches. In this context, the Landau- 
de Gennes [9j theory has been extensively used to predict properties of defects, but this 
theory is still mesoscopic in nature and makes no contact with particle interactions. An 
alternative is to use computer simulations, but these are generally time consuming for the 
study of defects. Therefore, the formulation of theories based on molecular approaches are 
needed. A microscopic theory, of the Maier-Saupe type, has been advanced [lOj, but it has 
some shortcomings; for example, it predicts all elastic constants to be equal, which causes 
different types of defects to have identical properties. This paper is devoted to exploring the 
consequences of another such theories, namely a simple version of density-functional theory 
(DFT) for hard anisotropic particles in two dimensions, which should give more realistic 
values for the size and energies of defect cores since the theory predicts different values for 
the elastic constants. 

DFT is ideal to study liquid-crystal defects, since it self-consistently gives the thermo- 
dynamic and microscopic structural properties of the inhomogeneous nematic fluid. One 
advantage of DFT over traditional approaches is that elastic constants, in particular the 
problematic surface elastic constants, and other phenomenological parameters, do not ap- 
pear explicitely in the theory, but only implicitely through interactions and distribution 
functions in a free-energy functional which is minimised (to all orders in the director spatial 
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derivatives). The defect-core structure appears naturally, and this is ideal since, contrary to 
the usual approximation within elastic theories of ignoring the defect cores in larger-sized 
nematic droplets, the contribution of defects cannot be ignored in nanocavities. One ques- 
tion is why details of defects should be important to understand large-scale configurations 
of the director field and defect motion. The microscopic approach enjoys some advantages 
whenever the relationship between bulk properties (such as elastic constants) and molecu- 
lar structure and interaction parameters is required. Knowledge of the detailed structure 
of defects will not be crucial to understand large-scale configurations and defect motion in 
stable nematics subject to boundaries or in nematic matrices where colloidal particles are 
embedded, but there are circumstances where this may not be so. For example, in the ki- 
netics of defect formation, re-organisation and anihilation, it may be important to know the 
core structure at short length scales. The microscopic approach can give useful estimates of 
free-energy changes, which are necessary to study the coarsening dynamics at a more mi- 
croscopic level using, for example, a relaxational dynamical equation. Also, the microscopic 
approach is essential at temperatures close to the clearing temperature, where defects act 
as nucleation seeds for the isotropic phase, and the structure and dynamics of the defects 
may be changing dramatically. 

In DFT the structure of the fluid is summarised by the local density and orientational 
distribution functions, which in turn may be used to obtain the more familiar nematic 
order parameter and local director field; these two are basic to describe nematic fluids 
containing defects in the director field. The connection between the two descriptions is done 
via the local one-particle distribution function, p(r,<f>), which gives the average number of 
particles at some position r with some orientation (on the two-dimensional plane). This 
quantity is obtained directly from DFT, and from it all interesting fields can be extracted, 
for example, the microscopic director field n(r), which is obtained locally as the direction 
where the orientational part of the one-particle distribution function presents a maximum 
(the macroscopic nematic director could be obtained by some coarse-grained average of 
the latter over some appropriate volume). Therefore the defect core region, along with the 
far neighbourhood of the singularity, can be analysed within a single framework based on 
particle interactions. The computational demands of the method are high, however, and in 
the present paper we restrict ourselves to the case of two-dimensional cavities of small radii 
(in the nm scale). 
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FIG. 1: Schematic of the two types of defects of topological charge k = +1 studied in this work, 
(a) Radial defect, with particles pointing on average along the radial direction, which excites splay 
distortion mode, (b) Tangential defect, with particles pointing on average along the tangential 
direction, which excites the bend distortion mode. Dark circular regions represent the defect cores. 
Lines are tangent to the director field. 



A defect is a singularity of the nematic director field n, characterised by a topological 
charge k, i.e. the number of turns of the director when the singularity is completely encircled 
6|. Elastic theory assumes smooth spatial variations of the director and therefore is not able 
to account for the structure of the singularity. In two dimensions the local elastic free-energy 
density can be written as 



1 2 1 2 

fei(r) = -h (V ■ h) + -k 3 \n x (V x n)\ 



(1) 



where ki,k% are elastic constants for splay and bend deformations (twist deformations are 
not possible in 2D). Let us consider the two defects with topological charge k = +1 depicted 
in Fig. (U which are called 'radial' (r) and 'tangential' (t). If <p is the polar angle of 
the position vector r, then the director field for the r defect is n = (cos 0, sin 0), and 
V-n = l/r, Vxn = 0, so that only splay deformations are involved. In the t defect 
we have n = (sin0, — cos0), V • n = 0, |n x (V x n)| = 1/r, and the only deformations 
involved are of bend type. A general deformation will involve both modes. Now, due to 
the singularity at the origin (location of defect), the elastic free energy within an arbitrary 
area containing the origin will diverge logarithmically: elastic theory fails here, and it is 
necessary to subtract this region by arbitrarily defining a core region, with free energy F n 
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and radius r n . The free energy within a circle or radius R will be: 

F r = jrA* tog A + 1*0 

F t = TrAfclogA + FW. (2) 

For 1? — > oo these energies diverge logarithmically, a situation that cannot arise in practice 
due to the presence of defects with opposite charge in the material. 

Little is known about the structure and properties of defect cores |7]. They are generally 
treated at a qualitative level, estimating the radius and defect-core energy in an approximate 
way [3]. Sometimes it is assumed these energies to be negligible compared with the elastic 
energy, and therefore defect cores are neglected altogether, a drastic simplification which can 
be severe if the system size is small. One of the first attempts to describe the defect core is 



due to Schopohl and Sluckin [111 ], who analysed a 1/2-disclination using Landau-de Gennes 
theory with the complete ordering tensor Q. This study demonstrated that the core of these 
defects does not consist of a region of isotropic material, but rather it is ordered along the 
disclination line (a possibility that does not exist in 2D). Later Monte Carlo simulations on 
a hard sphero cylinder model by Hudson y Larson |l2j] corroborated this prediction, and also 
found a new structure with a stable triangular nucleus for very elongated molecules. 

n 

In the only truly microscopic theory presented so far, Sigillo et al. [10[ used an approach 
based on a Maier-Saupe theory with an orientational distribution function, analysing discli- 
nation lines of charge k = +1 within a cylinder. This is a three-dimensional setup, while 
ours is a two-dimensional one. However, if one forgets about escape configurations, the 
director field should in this case exhibit the same kind of configurations as in our problem. 
The authors observed that the radius of the core decreases as the orientational order param- 
eter increases. Also, they examined radial and tangential defects and analysed their cores 
and their energies, obtaining that the two have the same size and energy. As the authors 
recognise this conclusion, which is certainly wrong, is due to the simplified interaction po- 
tential used, inherent in the Maier-Saupe theory, which predicts identical values for all the 
fluid elastic constants. 

Despite the reduced theoretical attention received, defect cores may play a very important 



role in many aspects of liquid-crystal science. Mottram et al. 13|, [14|] studied disclination 
lines of charge k = +1 and k = +1/2 near the isotropic-nematic transition in 3D, and 
explained the impossibility of heating a nematic material above a critical temperature T c > 



FIG. 2: Schematic of a hard discorectangle of total length L + D and width D (left), and cavity 
with a particle inside (right), showing definition of radial distance r, polar angle <f> and cavity radius 
R. Dashed particle has its centre of mass right at the cavity wall and cannot move further inside 
the wall (i.e. inside the cavity outer region). 

T IN (clearing point) which is below the limit of metastability of the nematic phase, due to 
the growth of the isotropic core. Defect cores properties may also be relevant in dynamical 
aspects such as defect motion Q, Q- 



In this paper we make a first attempt at calculating the properties of a defect core using a 
microscopic approach based on density-functional theory, using hard-particle interactions. 
One of our aims is to improve upon the results of Sigillo et al. [lO] by making more sensible 
predictions about the properties of the two types of defect cores investigated, namely with 
radial and tangential symmetries, in a circular cavity. In Sec. [Ill we briefly review the 
particle model and the DFT theory, together with the numerical approach and the bulk 
behaviour. Sec. II I II is devoted to the calculation of the elastic constants of the model. In 
contrast with elastic or Landau approaches, the DFT formalism does not require to specify 
which region is the core and which region is not the core, so that some criterion, similar to 
the Gibbs dividing surface in the statistical mechanics of interfaces, is needed to define the 
core. In order to analyse this problem it is necessary to compare the results from DFT with 
elastic theory, and this demands knowledge of the splay and bend elastic constants ki and 
&3. The values of these constants can be obtained within the same DFT framework. In Sec. 
II VI we present results for two types of point defects of charge k = +1 inside a circular cavity, 
placing emphasis on the size and energy of the defect core. Also, we discuss a configuration 
containing two k = +1/2 defects, for which no definite conclusions can be drawn (due to 
the cavities explored being too small) but some trends can be obtained. 
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II. THEORY 

In a previous paper [r| we have presented results for the structure, thermodynamics 
and structural (Frederiks) transitions of nematics confined into two-dimensional circular 
nanocavities using DFT. Here we use the same version of the theory in the same setup, but 
with an emphasis on defect core structure and energetics. Details of the theory were given 



in Ref. 16|; here we give a summary of the main features. 

The particle model used is the hard disco-rectangle (HDR), Fig. [2j which can be thought 
of as the projection of a sphero cylinder on a plane. A HDR particle has a rectangular section, 
of length L and a diameter D, and two semicircular caps at the two ends of the rectangle, also 
of diameter D. These particles interact via exclusion (i.e. configurations with overlapping 
particles are not allowed, but particles are not interacting otherwise), and can form a two- 



dimensional nematic at high volume fraction 17|. As interactions are hard, the temperature 



dependence is trivial, so the relevant intensive variable in the thermodynamics of this fluid 
will be the chemical potential (or, alternatively, the density). 

In DFT one writes an approximate free-energy functional F[p] in terms of the one-particle 
distribution function p(r,0), which can be split as p(r,(j>) = p(r)/(r,0), where f(r,4>) is 
the angular distribution function, and p{r) = J d(pp(r, 0) is the average local density. The 
free-energy functional is written as 

F\p]=F id \p]+F m \p]+F„ t \p] (3) 

with F id [p] the ideal contribution, 

(3FM = [ drp(r) {log [p(r)A 2 - l] - k^S^r)} , (4) 

J A 

where A is the total area of the cavity, A the thermal wavelength, and S Iot (r) the local 
rotational entropy density: 

S Iot (r) = -k r d<Pf(r, <P) log [27r/(r, 0)]. (5) 
Jo 

As usual, (3 = 1/kT, k being Boltzmann's constant. The excess part, F cxc [p], we write in 
terms of that of a reference fluid of locally parallel hard ellipses, which in turn is obtained 
exactly from that of a hard-disc fluid. HDR, ellipses and discs will be chosen to have the 
same particle area v and, in the case of HDR and ellipses, the same aspect ratio. These 
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conditions are sufficient to fix cry and a±, the diameters of the ellipses along the major and 
minor axes respectively, and from here a e , the hard-disc diameter, with of = cr\\<7±. In the 
following we will use x = L/D = 15 (which gives L = 3.346a e and D = 0.223a e ). The 
excess free-energy per particle of the hard-disc fluid is obtained from a theory due to Baus 
and Colot |ig|: 



f^Uv) = (ca + 1)^ + (P2 - 1) log (1 - V), (6) 

where C2 — 0.1280, rj = p$v is the packing (or volume) fraction, and po the mean number 
density. The excess free energy is then written as 

F exc [p}= dr d<f>p(r,<P)<p(r,<f>), (7) 

J A JO 

where the local free-energy per particle is 

x [ dr' [ d<t>'p{r', <t>')vUr ~ r', 0, <f>'). (8) 

J A JO 

w cxc is the overlap function of two HDR particles (equal to zero if particles overlap and unity 
otherwise). This expression is a variation of the Parsons-Lee 19) theory for homogeneous 
fluids of hard rods, or (from a different perspective) a variation of the Somoza-Tarazona [20I 
theory for inhomogeneous fluids of hard rods (both in three dimensions). Finally, F ext [p] is 
the contribution from the external potential (see Fig. [2]): 

F**\p]= dr rf0p(r,0)w cxt (r,0). (9) 

J A Jo 

The external potential acting on the particles will be chosen according to the type of favoured 
particle orientation at the cavity surface. In the case of the radial defect it is sufficient to 
use, as an external potential, a hard wall acting on the particle centres of mass: 

00, r > R, 

(10) 

0, r<R, 

where r is the radial distance measured from the centre of the cavity, and R is the radius of 
the circular cavity. This choice is known to favour homeotropic (i.e. perpendicular to the 
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wall) orientation of the fluid director next to the wall [if]], thus inducing an r-type defect. 
For the tangential defect a different choice is necessary (see Sec. IIV 2j) . 
The angular distribution function f(r,(p) is parameter ised according to 

a(r) cos 2[rf>— «/>(r)] 

f{rA) = —« , (11) 

d4>e a{r)cos2 ^ 

where the field ip(r) is the local tilt angle of the nematic director, measured with respect to 
the x axis, and a(r) is a variational function, related with the local nematic order parameter 
q(r) by 

q(r)= [ rf0/(r,0)cos{2[0-^(r)]}. (12) 
Jo 

These equations allow one to describe the configuration of the fluid by means of the three 
local fields {p(r), q(r), ip(r)}; in the following we will use the local packing fraction, rj(r) = 
7TO"gp(r)/4, instead of the local density p(r), as basic density variable. Finally, for a cavity 
of fixed radius R, we impose on the system a constant chemical potential p, and minimise 
the cavity grand potential 

h[p] = F[p] -pi dr ^ d<p P (r, <P) (13) 

J A Jo 

with respect to variations of the variables defined above. To obtain the minimum the two 
dimensional space xy is discretised into a square lattice with spacing Ax = Ay = 0.089<7 e , 
with mesh points (xi,yj), representing 40 points in a particle length L + D. The circular 
surface is approximated by a zigzag line. The trapezoidal rule was used to calculate spatial 
integrations, while angular integrals were approximated using Gaussian quadrature with 
30-40 roots. The field variables {rj(r), q(r), ^(r)} were discretised as i]^, qij and and 
the free-energy functional was minimised using the conjugate-gradient method. 

This model presents a bulk isotropic-nematic phase transition for packing fraction ^ IN = 
0.257 and reduced pressure pv /kT = 0.98 (estimates from simulation [l7] give 7] m = 0.363). 
The transition is of the second order. 



III. ELASTIC CONSTANTS 



In order to compare with elastic theory, we need some criterion to define the boundary 
of the defect core. In this respect it will be useful to compare the free-energy densities from 
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DFT and elastic theory since, from this comparison, we can locate the boundary separating 
defect core from the outside region (where elastic theory should be valid) as the distance 
where both densities coincide. We will see that this definition is somewhat arbitrary, as 
it relies on our definition of how close, numerically speaking, the two free-energy densities 
should be. We come to this point later. For the moment, we note that the elastic free-energy 
density contains elastic constants that we must be known in advance. These constants have 
to be calculated within the same DFT scheme: the DFT free-energy density will smoothly 
tend to the value predicted by elastic theory provided we use the values for elastic constants 
predicted by DFT. Then we calculate separately the elastic constants k\ and ks in the 
framework of DFT. The expressions for the elastic constants are: 

4? 7o Jo Jo 

h = -^^ F 64 F d^pw^vu^n 

4r 7o Jo Jo 

(14) 

where r] = p v is the bulk packing fraction. These constants are evaluated at the uniform 
nematic (no spatial inhomogeneities or director distortions). The one-particle distribution 
function is then p(r,0) = p(0) = p f '(</>). In the expressions above, p'(4>) is the derivative 
of the one-particle distribution function with respect to the tilt angle, p'(0) = dp/dip = 
Podf/dcj), and where we defined 

Vij((j), 0') = / drv cxc (r, 0, <j}')xiXj. (15) 

J excl. area 

The area integral over r is extended over the area of exclusion of two particles. Details on 
how these expressions are obtained are to be found in the Appendix. An alternative and 
equivalent way to obtain the elastic constants is to use the same confinement setup (circular 
cavity) defined above and impose a given director field with pure splay or bend deformations 
(Fig. [[]), setting the density and nematic order parameters to the corresponding bulk values. 
Now if the free-energy density, as given by evaluation of the functional, is represented 
along any one of the cavity diameters (there is azimuthal symmetry), we can extract the 
elastic constants by comparing with the radial dependence predicted by elastic theory in 
the intermediate region (far from both the cavity centre and the cavity surface). We have 
seen already that the radial dependence is ~ 1/r 2 in both cases. Since no minimisation is 
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Afi/kT rjo q fa/kT k 3 /kT 

0.15 0.270 0.27 0.08 0.11 

0.75 0.303 0.63 0.47 1.13 

1.75 0.360 0.83 0.93 4.05 

2.75 0.410 0.90 1.27 8.72 

4.25 0.470 0.94 1.70 17.8 

5.25 0.503 0.96 2.01 25.2 

6.25 0.533 0.97 2.38 33.8 



TABLE I: Bulk properties of nematic fluid of HDR of aspect ratio x = 15, as obtained from DFT. 
Afi/kT is the excess chemical potential with respect to the isotropic-nematic coexistence value, in 
units of thermal energy kT; t]q and q are the packing fraction and nematic order parameter; and 
ki/kT, k^/kT are values of elastic constants, also in units of thermal energy. 

implicit in this method, one can use very large cavities [R ~ 100(L + D)] so that the elastic 
constants can be obtained with accuracy. 

In Table Q] and Fig. [3] we provide values for ki and k 3 (as expected, the two strategies 
to obtain the elastic constants explained above give the same results, except for some tiny 
differences that come from the numerical accuracy of angular and spatial integrals). The 
values of the elastic constants are zero at the bulk transition. k 3 is always larger than hi, 
and their difference increases with density: when rjo ~ 0.4 the difference is almost an order 
of magnitude, which means that bend deformations are more costly energetically than splay 
deformations. This is an important point, as most studies based on elastic theory assume the 
one-constant approximation k± = k 3 . In our case (HDR particles with aspect ratio x = 15) 
this approximation ceases to be valid even very close to the isotropic-nematic transition [see 
Fig. [31(d)] . In confined nematics under strong geometric restrictions such as the one studied 
here, nematic order is very frustrated and stable nematic configurations are only obtained 
for conditions deep into the bulk nematic stability region (i.e. and considerably far from the 
bulk transition); this means that the one-constant approximation will be very inaccurate. 
For particles with lower aspect ratios this problem will become less acute. 

One consequence of this problem can be seen in the paper by Bates [21[ , where the nematic 



12 




0.2 0.4 0.6 0. 

q 
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FIG. 3: Elastic constants of a nematic fluid of HDR particles with aspect ratio x = 15, as obtained 
from DFT. (a) Splay elastic constant k\ in thermal energy units kT as a function of packing fraction 
T). (b) Bend elastic constant k^ in thermal energy units kT as a function of packing fraction, (c) 
Splay (continuous line) and bend (dashed line) elastic constants as a function of nematic order 
parameter q. (d) Ratio of elastic constants as a function of packing fraction. 



ordering of hard sphero cylinders lying on the surface of a sphere is examined via Monte Carlo 
simulation. Geometry forces the creation of four defects of charge +1/2. However, analysis 
based on the one-constant approximation predict that the defects are located at the vertices 
of a tetrahedron, while the simulations show that they are in fact distributed along a great 
circle: in this way the director field arranges itself in a way such that splay distortions are 
maximised, while bend distortions, much more costly energetically, are minimised. 

Another observation of our calculations concerns the bulk isotropic-nematic transition. 



This transition has been studied by Bates and Frenkel [17] b y M onte Carlo simulation. 
Assuming the transition to be of the Kosterlitz-Thouless type [22j, and also that the two 
elastic constants are equal, the transition should occur when the elastic constant reaches the 
critical value k c = 8kT/ir. Using our values for the elastic constants and taking the average 
k — (ki + k 3 )/2, we obtain i] m = 0.36, in perfect agreement with the simulations. 
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IV. RESULTS 



In this section we analyse various types of defects. We start with the radial configuration, 
r, where only splay-type director distortions are present and there is a central defect of 
topological charge k = +1. There follows the case of charge k = +1 but with a tangential, 
t, director field. Finally, we will consider point defects with charge k = +1/2. 



1. Radial defect with charge k = +1 
In Ref. 



161 ] we found that this defect can only be stabilised at low chemical potential, 
close to the bulk isotropic-nematic transition. As the chemical potential is increased, the r 
configuration becomes metastable and the central k = +1 defect splits into two k = +1/2 
defects. However, it is possible to impose the r configuration by preparing the system so that 
the director is forced to always point radially, keeping the director field unchanged during 
the conjugate-gradient minimisation. 

An example is given in Fig. HI where the packing-fraction profile along one diameter is 
displayed. The different cases shown correspond to increasing cavity radius, from R/(L + 
D) = 1.98 to 7.98. Calculations are presented for fixed relative chemical potential A^/kT = 
2.75, where A/z is referred to the value of fi at the bulk isotropic-nematic transition. Three 
well-defined regions can be seen. In the central region a marked depletion in number of 
particles is observed, which corresponds to the defect core. In the neighbourhood of this 
region the density is quite constant, and as the inner surface of the cavity is approached 
a local minimum appears, followed by a sharp density increase due to surface adsorption. 
Since we are mostly interested in the defect core, to minimise the effects of the surface on the 
core properties we need to consider as large a cavity as possible. Our present computational 
capabilities limit the radius of the cavity to R ~ 10(L + D). However, a simple inspection 
of the profiles seems to indicate that the surface effects are relatively weak, even in small 
cavities. This figure clearly demonstrates that the size of the defect core is well defined even 
for cavities of small radius [say R > 3.2(L + D)]. 

In the neighbourhood of the defect core there is a region dominated by elastic effects. If 
the defect were very far from any surface this region would extend up to the surface, but the 
question is: is it possible to obtain a truly elastic regime for small cavities such as the ones 
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FIG. 4: For a radial defect of charge k = +1, local packing fraction i) as a function of radial 
distance r from the cavity centre (in units of particle length L + D), along an arbitrary cavity 
diameter, at relative chemical potential Afi/kT = 2.75. Lines correspond to different cavity radii: 
R/(L + D) = 1.98,3.17,3.98,4.98,5.98 and 7.98 (continuous and dashed lines alternate for a better 
visualisation) . 

investigated here? To answer this question, we focus on the (grand-potential) free-energy 
density inside the cavity, uj(r), defined by 



In Fig. [5] the free-energy density is plotted as a function of radial distance from the centre 
of the cavity, for a cavity radius R = 7.98(L + D). The elastic free-energy density / e i(r) = 
k\/2r 2 is also included; to obtain this energy, the value for the k± elastic constant was 
taken from the DFT calculations (Table [I]). Of course both free energies disagree in the 
central region of the cavity (where the free-energy density from elastic theory diverges at 
the singularity) and in the region close to the surface. However, there is an intermediate 
region, in the interval r e <r <r s , [withr e ~ 2.0 — 3.0(L+D) and r s ~ 5.0— 5.5(L+D)] where 
the agreement is quite good; this is a signature of the elastic region. The conclusion that an 
elastic regime can indeed be defined was also reached by Sigillo et el. [h]] in their Maier- 
Saupe approach and indirectly in Landau-de Gennes approaches |ll|. The free-energy 
density may be used to loosely define a defect-core size in terms of the radial distance at 
which the free-energy density begins to behave as ~ r~ 2 (of course this is an ambiguous 
definition that, in practical terms, does not affect the numerical values of the defect-core 
properties significantly). In the case of Fig. [5] we obtain a size 2r e ~ 4 — 5(L + D); this 




(16) 
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1 2 3 4 5 6 7 
r / (L + D) 

FIG. 5: Grand-potential density u, in units of a e and kT, as a function of radial distance (in units 
of particle length L + D), along an arbitrary diameter of a cavity with R = 7.98(L + D) (continuous 
line). The relative chemical potential is Afi/kT = 2.75. The dashed line is the elastic free-energy 
density according to elastic theory. r e and r s are approximate radial distances for the boundaries 
of the elastic region. The inset is a zoom of the central region. 



should be a few times the correlation length £, which is in agreement with calculations on 



three-dimensional defects by Landau-de Gennes theory [11]. 

The properties of a cavity of radius R = 7.98(L + D) are summarised in Fig. [6] The 
profiles of the nematic order parameter, Figs. [6](a) and (b), indicate that the core radius 
decreases as the chemical potential fi increases. To quantify this effect more precisely and 
analyse the depletion of the order parameter that occurs inside the central region, we have 
denned two additional measures of the defect-core radius, and rn \ as the inflection 
points in the nematic order parameter and the density profiles, respectively; these two 
quantities do not coincide with, but should behave like, the energy-based measure r e as 
thermodynamic conditions are varied). In Fig. M^c) we plot these quantities as symbols. 
Both have a similar behaviour: they decrease quickly with \i and saturate at high chemical 
potential, with the inflection point of the nematic order parameter saturating a bit earlier. 

The fact that the core radius decreases with \i does not mean that its effects propagate 
to a smaller region; in fact, the result is quite the opposite. The difference between bulk and 
core densities increases with /z, while the intermediate region extends to larger distances. 
However, the effective core radius r e (where the free-energy density differs significantly from 
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FIG. 6: Some properties of a nanocavity of radius R = 7.98(L + D). (a) Local packing fraction rj as 
a function of radial distance r from the cavity centre (in units of particle length L + D) for various 
values of relative chemical potential: A^i/kT = 0.75, 1.75, 2.75, 4.25, 5.25 and 6.25 (from bottom to 
top) and for a radial defect with k = +1. (b) Nematic order parameter q as a function of radial 
distance r from the cavity centre (in units of particle length). Profiles as in panel (a), (c) Core 
radius r n in units of particle length as a function of relative chemical potential. Filled circles: . 
Open circles: r^ 2 \ (d) Core energy F£' as function of splay elastic constant k\, both in thermal 
energy units. Error bars were calculated with the two choices r e = 2.5(L + D) and 3.0(L + D) for 
the upper limit in the integral of Eqn. (|17p . The straight line is a linear fit. 



the elastic one) is usually in the interval 2.5 — 3(L + D), largely independent of /i for large 
A/i. 

The successful identification of an asymptotic elastic region and the ensuing possibility 
of defining a defect-core boundary allows us to associate a free energy F n with the defect 
core. This we do by integrating the excess of grand-potential density over a uniform fluid 
at the same chemical potential inside a circle of radius r e (where elastic behaviour sets in): 

F n = 2n / drruj(r). (17) 
Jo 

This energy is represented in Fig. [6]^d) as a function of the elastic constant k\. The calcu- 
lation has been done using r = 2.75(L + D) as a cut-off distance, but calculations were also 
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done using 2.5(L + D) and 3.0(L + D) to see the effect of changing the cut-off; error bars 
in the data correspond to these two limits. As can be seen in the figure, the differences are 
very small. 

In phenomenological treatments it is usual to assume that the free energy of a disclination 
core of charge k is F n = k 2 irk, where k is an elastic constant [3] in the one-constant 
approximation k = k\ = k^) which, for the radial defect, is k\. Our DFT results give 
support to the linear relation between Fn and fci, but the slope (obtained by a linear fit) is 
equal to 12.2, which is four times larger than that predicted by the phenomenological theory 

(r) 

for a point defect of charge k = +1. The dependence of Fn on chemical potential is also 
linear, with a slope of 4.09 (not shown). 



2. Tangential defect with charge k = +1 

In this case the director field only supports bend distortions, as shown in Fig. [TJ To 
stabilise such a structure we need a surface potential that favours planar anchoring, i.e. 
particle orientations tangential to the surface. We use the following model for external 
potential: 



oo, r > R 



(18) 

V r o cos2(0-^)e- Q( - R - r) , r<R, 



where Vq is the surface strength. For Vq large enough, the surface favours tangential anchor- 
ing; we have checked that this is the case, e.g. for Vo = 0.7kT and a = 1.08(L + £)) _1 . These 
are the values we will be using in the following to study a defect with tangential anchoring. 

The inclusion of an external field with an exponential decay means that the surface 
interacts with the fluid at longer distances than in the previous case. An additional feature 
is that, since the splay elastic constant is smaller than the bend elastic constant, the size of 
the defect core is larger. Both these effects play against the possibility of reaching the elastic 
regime in the region between the defect and the surface. Therefore, much larger cavities are 
needed. Our computational limit is R ~ 15(L + D), which is not large enough to obtain 
reasonably accurate estimates of the core energy, for example. The only safe conclusion is 
that this energy is significantly larger than that of the radial defect. 
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FIG. 7: (a) Local packing fraction rj as a function of radial distance r from the cavity centre 
(in units of particle length L + D) for various values of relative chemical potential: Afx/kT = 
0.75, 1.75,2.75,4.25,5.25 and 6.25 (from bottom to top) and for a tangential defect with k = +1. 
(b) Nematic order parameter q as a function of radial distance r from the cavity centre (in units 
of particle length). Profiles as in panel (a). All data pertain to the case R = 7.98(L + D). 

Despite this problem, it is instructive to study the structure of the defect core in a 
qualitative way. Fig. [7] shows the order parameter profiles in the defect-core region for a 
cavity of radius R = 7.98(L + D) (for larger cavities the profiles will be slightly different). 
The nematic order parameter behaves similarly as in the previous case, save the different 
size. The density has a pronounced maximum at the core centre. The size of the central 
depleted region is larger than one particle length, allowing for a higher particle concentration 
inside the defect core. 

3. Defect with charge k = +1/2 

The present geometry can also be used to explore a more interesting defect with 

charge k = +1/2. The minimum-energy state contains two defects of charge k = +1/2, 
separated by a distance d , when the cavity radius is sufficiently large, but the analysis is 
more complicated here, as two additional minimisations are required: a partial one with 
respect to the 'fast' variables at fixed defect separation, and a minimisation with respect 
to the defect separation d,Q, which is a slow variable. As a result, the computation time 
increases by an order of magnitude. The practical consequence is that the maximum radius 
of the cavity that can be analysed is reduced, and the task of splitting contributions of defect 
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FIG. 8: Contour plots with respect to xy coordinates for local packing fraction r] (upper row), 
nematic order parameter q (middle) and director tilt angle Vl/ (lower) for two configurations with 
relative chemical potentials = 1.75/cT (left column) and Afj,/kT = 6.75 (right column) in a 
cavity of radius R = 3.18(L + D). xy coordinates in units of D, and tilt angle is given in degrees. 



cores from the rest becomes harder. 

In Fig. [8] we have plotted the order parameters for a configuration with two k — +1/2 
defects, at two different chemical potentials. The configurations were obtained by minimising 
the functional in a cavity of radius R = 3.18(L+D) (only a region of size 30 x 30D 2 containing 
the two defect cores is shown). In the left column local packing fraction rj (top), nematic 
order parameter (middle) and tilt angle (bottom) are shown for the case A/i = 1.75kT. 
The right column shows the same profiles when the chemical potential is increased to A/i = 
6.75&T. In the corners of the density plots the structure has radial symmetry: this is a 
surface effect. This is an indication that larger cavities may be necessary for a more detailed 
study. We can clearly see that the core size decreases considerably as the chemical potential 
is increased (this effect is more visible in the nematic order parameter). Another remarkable 
effect is the loss of radial symmetry of the defect as the chemical potential is increased. The 
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FIG. 9: (a) Radius r n as a function of relative chemical potential A/x, and (b) core energy F n as a 
function of elastic constant k\ for a defect of charge k = +1/2. 



profiles in the left column (low chemical potential) have an almost radial symmetry with 
respect to the defect core (save the tilt angle, obviously). As /x is increased (right column), 
the core shrinks in all directions, especially along the direction joining the two defects, where 
the director field is constant. 
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Unfortunately, we have not been able yet to study cavities large enough for the core 
structure and the surface structure to relax completely to the elastic limit, but qualitative 
estimates can be obtained for the relevant properties of the core. This can be seen in Fig. 
In panel (a) we plot the average radius of a defect core as a function of the chemical potential. 
Similar to the radial defect, the average radius r n has been defined as the inflection point of 
the density profile, averaged over all directions (since here there is no angular symmetry). We 
can see that there is a rapid decay as /i increases, and r n levels off at a value approximately 
equal to half the value for the radial defect [see Fig. HE(c)]. Therefore, the core area, which 
is proportional to r^, is about four times less in the case k = +1/2 than in the case k — 1. 
In panel (b) of the same figure, we plot the defect core energy as a function of the elastic 
constant k\. The straight line is a linear regression with slope m = 5.04, approximately a 
factor 2.5 smaller than in the k = +1 defect. We could expect a factor 2 beforehand, since we 
know in advance that, even in very small cavities, the structure with two k = +1/2 defects 
is more stable than that with a single, radial defect. In this case the largest contribution to 
the free energy comes from the defect cores, and therefore the energy of both cores plus their 
repulsive energy must be at most equal to the energy of the radial defect. We mentioned 
before that the defect energy is generally taken to be F n = k 2 Trk and we would expect a 
factor 1/4 in the case k = +1/2 with respect to the radial defect (the size is approximately 
four times smaller). The behaviour of the director field in the core region exhibits bend-like 
distortions (V x n ^ 0) for k = +1/2, which do not appear in the radial defect and can be 
the origin of such a difference. 

V. SUMMARY AND CONCLUSIONS 

In summary, we have studied the core properties of a defect in 2D, using a DFT (micro- 
scopic) model, free from fitting parameters, that includes consistently variations in density 
and nematic order parameter. All cases studied predict the formation of an isotropic region 
in the defect, something to be expected in 2D. The core free energy is proportional to the 
elastic constant, with a varying proportionality constant that depends on the type of defect 
studied; this is due to the different energetic cost associated with deformations of splay and 
bend type. The size of the defect cores is on the order of a few particle lengths, and decreases 
as the nematic ordering of the surrounding fluid increases. The core size saturates for strong 
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nematic ordering. 

The study of a 2D defected nematic fluid, used here mainly for computational reasons, 
may be useful to understand 3D phenomena in the physics of defects. As mentioned in 
the introduction, knowledge on the structure and energetics of defects may be important in 
dynamical problems, such as defect formation or nucleation and coarsing of the nematic and 
isotropic phases. This structure may be changing in time and a microscopic approach may 
be helpful to follow the dynamics via relaxation equations that involve the gradient of a free 
energy. A natural extension of our work therefore involves the study of 3D nematic fluids 
and their defects. Schemes where the microscopic and mesoscopic approaches are combined 
would also be useful in the above-mentioned problems, and work along this avenue is under 
way in our group. 
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APPENDIX A: CALCULATION OF ELASTIC CONSTANTS 



astic constants of a 3D nematic liquid crystal were derived by 



Expressions for the e 

Poniewierski and Stecki 23j using a direct correlation function route. ^From these expres- 
sions it is easy to write the corresponding 2D expressions. Obtaining the direct correlation 
function of the theory (which, in our Onsager-type theory, is basically the Mayer function) 
one can obtain explicit expressions in terms of integrals over the excluded area and the 
orientational distribution functions. Here we present an alternative derivation, valid only in 
2D, in terms of expansions in the local tilt angle ip{r). We start from Eqns. (J7|) and (JS]) for 
the excess free energy of a nematic with constant density, expressed explicitely in terms of 
the tilt angle: 

F CXC [ P ] = jdrj dMd>-Tp(r)) 

x f [ dr'd<j>'v m (r - r', 0, <f>')p[<f>' - V(r')] (Al) 
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Now we expand the second local density p(4>' — ip(r')) in ip(r') around 
Ar — r' — r, we have: 

il){r') = if)(r) + Ar ■ V r ^ + - [Ar • V r ] 2 ip H 



Then we expand the density: 



'-V(r')] = P W - ^(r)] + ^ 



(rO-^(r)] 



+ 



1 <9 2 p 



2 



(r') - ^(r)] 2 + 



Substituting (1A2|) and keeping terms up to third order in the gradient: 

dp 



p[ (f) '-^(r l )}=p(<p'-^(r)) + 



dip 



x ( Ar • V r ip + - (Ar • V r y ip + O (V r ^T 



+ 



1 d 2 p 



2 a^ 2 



Ar • V r ip + - (Ar ■ V r ) 2 ^ + O (V r ^) 3 



+ ... 



p[0' - V(r)] + p;[0' - V(r)]Ar • V r ^ + 



(Ar- V r ) 2 ^ 



+- 



^(r)] (Ar • V r ^) 2 + 
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The first term gives the free energy of the undistorted nematic (since ip(r) is supposed to 
be a slowly varying field). The elastic free energy is then: 

x drd(j)p((f) - ip(r)) / / dr'd(f)'v cxc (r' ,(f),(f)') 



x - ^r)]r' • V r ^ + l -p^W - ^(r)] [r' • V r ] 2 ^ 

+ ^;(0'-^))[r'.V^] 2 .} (A5) 

We can take ip(r) = in the argument of the density profiles and its derivatives; the elastic 
free-energy density is then: 

/d(r) = J d<t>p{<t>) j j dr'd<t>'vUr'AA') 

x \p'^') (r • V r ^) + ip;(0') (r' • V r ) 2 V 

+^;(0')K-V^) 2 } (A6) 

Now we note that 

y drVc(T-^»' = 0, (A7) 

due to the symmetry V(r, 0, 0') = V(— r, <p, 0'), and the term linear in the gradient of ip{r) 
vanishes, as it should be. Defining the dyadic 

V(0,0') = J drv m (r,<i>,<j>')rr, (A8) 
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we get 

*exc(?7o) 



x E// Wt'VfrM) J drp^-^r)) 
x { P ;[0' - V(r)]9^V(r) + ~ i>(r)}d^(r)d^(r)} 

Now we integrate the term with the second derivatives by parts: 

J drp[(j>-ip{r)]p'^[(j>' - ij)(r)]dpyip(r) 



dr{p;[0-^(r)] P ;[0'-^(r)] 



,3=const. 



+ P [0 - i>{r)]p'^>' - ij(r)}}d^(r)d^(r). 
The pp" terms in (1A9|) cancel out, and the surface term is neglected. Then: 

^exc(??o) 



Pi 



P ;(0 - ^(r)) P ;(0' - ^(r))^(r)^(r), 
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which is the searched-for expression. In terms of the unit vector Cj along the particle axis, 
we have Co ■ n = cos (0 — -0), and 



and therefore 



P^O - = P'{& ■ n) sin (0 - -0), 



P ;(0'-^)=p / (^-n) sin (0' 



dff (a) • n) = uj a d/3n a = sin (0 - ifidpip, 



d 7 (uV • n) = uj' s djns = sin (0' — ip)d.yip, 



P^tt - i>)9/Pp = P'(& ■ n) sin (0 - ip)d p ip 



p'{Cj ■ n)} j uj a dpn a , 

a 



P'^' ~ ^)9^ = p'{&' ■ n) sin (0' - tfj)d y ip 



so that in terms of the gradients of the nematic director: 

*exc(%) 



a/375 



/ / d(j>d<j>'V l 



x p'(a; ■ a)p'(u)' ■ h)L0 a u' s dpn a d 7 ns. 



Since the direct correlation function of our model is 

*exc(?7o) 



c(r, U3, uj 
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(A12) 



(A13) 



(A14) 



our expressions for k\ , k 3 coincide with the general ones by Poniewierski and Stecki 23[ using 
a direct correlation function route. 
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Now in 2D the Frank elastic free energy contains only splay and bend distortions 

f cl (r) = h 1 (V-h) 2 + ^k 3 \Vxn\ 2 
Splay and bend are given, respectively, by the deformations: 

<*•*>•- (£ + S0* 

= (d x n x ) 2 + (d y n y ) 2 + 2 (d x n x ) (d y n y ) , 

, x _ f dn y dn x \ 2 
\dx dy J 

= (d x n y ) 2 + (d y n x ) 2 - 2 (d x n y ) (d y n x ) , 

and therefore: 

h = / / 

Here it is assumed that the director goes along the x axis. 
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